Electronic polarization in pentacene crystals and thin films 



m 
o 
o 

(N 



(N 



> 

in 
in 
m 
O 
m 
O 

c3 



O 

o 



13 



E.V. Tsiper 1 and Z.G. Soos 2 
1 Department of Chemistry, Rutgers University, 610 Taylor Rd., Piscataway, NJ 08854 
2 Department of Chemistry, Princeton University, Princeton, NJ 08544 

Electronic polarization is evaluated in pentacene crystals and in thin films on a metallic substrate 
using a self-consistent method for computing charge redistribution in non-overlapping molecules. 
The optical dielectric constant and its principal axes are reported for a neutral crystal. The polar- 
ization energies P+ and P- of a cation and anion at infinite separation are found for both molecules 
in the crystal's unit cell in the bulk, at the surface, and at the organic-metal interface of a film of 
N molecular layers. We find that a single pentacene layer with herring-bone packing provides a 
screening environment approaching the bulk. The polarization contribution to the transport gap 
P = P + + P_, which is 2.01 eV in the bulk, decreases and increases by only ~ 10% at surfaces 
and interfaces, respectively. We also compute the polarization energy of charge-transfer (CT) states 
with fixed separation between anion and cation, and compare to electroabsorption data and to sub- 
molecular calculations. Electronic polarization of ~ 1 eV per charge has a major role for transport 
in organic molecular systems with limited overlap. 
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I. INTRODUCTION 

Prospective molecular electronic devices are based on 
organic films deposited on metal or semiconductor sur- 
faces [1-3]. Due to limited mobility, organic devices are 
typically restricted to thin (10 — 100 nm) films. Varia- 
tions in electronic polarization energies of charge carriers 
near surfaces and interfaces affect the transport states [4] 
and hence device function. Wc have recently addressed 
electronic polarization at surfaces and in thin films [5] 
using a self-consistent approach in the limit of vanish- 
ing intermolecular overlap [6], when each molecule ex- 
periences the non-uniform electrostatic potential of all 
other molecules. For the prototypical hole conductor 
[3], perylenetetracarboxylic acid dianhydride (PTCDA), 
the calculated transport gap varies by 500 meV between 
monolayers and thick films, in agreement with photoclcc- 
tron and tunneling electron spectroscopy [5]. PTCDA 
showed that electronic polarization is accessible to self- 
consistent computation in crystalline thin films. 

In contrast to PTCDA, whose molecules form one- 
dimensional stacks and lie almost flat on metallic sub- 
strates, many organic molecular crystals of interest have 
a herring-bone packing, with molecules oriented across 
molecular layers as sketched in Fig. 1. Both thioplicncs 
[1] and acenes [2] have herring-bone packing and are suit- 
able for thin film transistors. The molecules in crystalline 
thin films are almost upright and charge transport is pref- 
erentially parallel to the surface. 

In this paper, we analyze electronic polarization in pen- 
tacene, considering it as a representative of a wider class 
of organic materials with electronic applications. We cal- 
culate electronic polarization in the bulk to obtain the 
optical dielectric constant, the transport gap for gener- 
ating an electron-hole pair at infinite separation, and the 



electrostatic binding between molecular ions at fixed sep- 
aration in the lattice. We then consider the experimental 
situation in Fig. 1 of thin films on a metallic substrate 
to compute electronic polarization at the surface and at 
the organic-metal interface, assuming idealized films with 
structure identical to the bulk. Crystalline thin films 
exhibit multiple phases that, except in monolayers, are 
close to the bulk. 




Fig. 1 Schematic molecular packing in a pentacene layer 
on metal. Image charges below the metal surface are also 
sketched. 

The long axis of pentacene or a-sexithiophene is almost 
normal to the surface. Image charges across the interface 
in Fig. 1 then act in the direction of greatest polarizabil- 
ity, which results in contrasting electronic polarization in 
pentacene and PTCDA films. Electric fields normal to 
the surface induce large dipoles in pentacenes, but these 
induced dipoles are parallel and hence repulsive. There is 
competition between charge redistribution due to image 
charges and induced dipoles. Our self-consistent calcula- 
tions of charge redistribution indicate that even a single 
molecular layer provides a screening environment that re- 
duces repulsion between induced dipoles. The improved 
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iterative procedures of Section II are necessitated by the 
competition between image charges and induced dipoles. 

The paper is organized as follows. In Section II we 
consider the system of linear equations for electronic po- 
larization and introduce alternative methods based on a 
minimum principle when the efficient iteration scheme 
of [6] fails. Section III contains results for bulk crys- 
tals, including the optical dielectric constants, polariza- 
tion energies of cations and anions, and the binding en- 
ergies of ion pairs. In Section IV we consider pentacene 
films of 1 — 10 molecular layers and report polarization 
energies at the surface and at the metal-organic inter- 
face. The self-consistent pentacene potential at the sur- 
face illustrates the competition between image charges 
and induced dipoles. Section V contains discussion and 
conclusions. 



II. SOLUTION OF LINEAR EQUATIONS 

The zero-overlap approximation reduces electronic po- 
larization to charge redistribution on molecules in the 
electrostatic potential of the solid. The crystal structure 
fixes all distances. We neglect lattice relaxation whose 
polarization is estimated to be an order of magnitude 
smaller in organic molecular crystals [7,8]. The crystal 
potential at atoms is readily introduced as site energies 
in semiempirical theory, which then provides a practi- 
cal approach to computing self-consistent atomic charges 
and induced dipoles [6]. We have eight linear equations 
per atom and need large (10 3 ) clusters of molecules for 
long-range Coulomb interactions. The large number of 
equations precludes the use of algorithms based on trans- 
formations of dense matrices. For example, a cluster of 
2000 PTCDA molecules, as in [6], leads to 608,000 lin- 
ear equations for 38 partial atomic charges and 3 x 38 
components of induced atomic dipoles for each atom of 
each molecule, and the same number of potentials and 
potential gradients. 

We review briefly some approaches for solving linear 
equations that appear in polarization problems. We first 
recast these equations in a formal but transparent form. 
Let q be the vector describing charge redistribution in a 
cluster of interest. Its components may include partial 
atomic charges for each atom of each molecule, as well 
as induced atomic dipoles. Higher atomic multipoles can 
be included in q as well. Let p(°) be the "unrelaxed" 
charge distribution, that is the charge distribution in the 
individual molecules in the gas phase, the actual charge 
distribution thus being p^ + q. 

For zero overlap, charge conservation for individual 
molecules leads to simple constraints on the components 
of q, meaning that q belongs to a certain subspace Q of 
interest. We note that p^ does not necessarily belong 
to Q since the source charges on ions do not sum to zero. 



The polarization problem is defined in terms of the 
"state vector" q as minimization of the total energy E(q), 
which can be formally written as 

E=\{qGq)+ l -(p^+q\V\p^+q). (1) 

We use bra- and ket- notation to denote matrix vector 
products. The first term is the energy of non-interacting 
molecules. A symmetric positive-definite stiffness matrix 
G ensures that q = in the absence of interactions. It 
results from solving Schrodinger's equation for a single 
molecule in non-uniform external field. The stiffness ma- 
trix is block-diagonal, with separate blocks corresponding 
to individual molecules or ions. 

The second term describes intcrmolecular interactions. 
The interaction matrix V contains charge-charge, charge- 
dipolc, and dipolc-dipole interactions, and can also be 
made symmetric. The matrix V is indefinite, reflecting 
that Coulomb systems are always unstable. In contrast 
to the stiffness matrix, the diagonal blocks of V are zero, 
since no molecule interacts with itself. 

Expression (1) shows a potential for instability, since 
it may have no minimum when V is sufficiently large. 
In the organic molecular crystals that we consider, V is 
small enough and only results in a shift of the minimum 
of E(q) from q = 0. 

Differentiating (1) with respect to q we find a set of 
linear equations on q: 

Gq + Vq = -Vp (Q \ (2) 

subject to the constraint that q belongs to Q. While this 
form of polarization equations is simple conceptually, an- 
other form is often more practical. By introducing a set 
of generalized potentials p = — Gq we recast (1) in a dual 
form: 

p=V(/V+q) (3a) 
q = -Up, (3b) 

where II = G _1 . These are formally written Eqs. (8) 
and (15) of [6]. The physical meaning of Eq. (3) is sim- 
ple. Instead of a one-step minimization of E(q) as in 
(2), we optimize the charge distribution of each molecule 
individually in the self-consistent external field p of all 
other molecules. 

The components of the vector p contain potentials 
and potential gradients (and possibly higher derivatives) , 
that couple to atomic partial charges, atomic induced 
dipoles, and higher atomic multipoles. The symmet- 
ric positive-definite polarization matrix II has the same 
block-diagonal structure as G and describes the polar- 
ization response of individual molecules to the external 
field, given by the components of p. The inverse of G is 
well defined in the subspace of interest Q. 
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A. Iterative solution of polarization equations 

The dual form of Eqs. (3) suggests a natural iterative 
approach. We start with the bare charge distribution 
p(°) and find "bare" potentials p^ = V p^> . We then 
repeat updating charge distribution using Eq. (3b) and 
recalculating p using Eq. (3a) until convergence. 

When convergent, this procedure is usually fast. How- 
ever, it is often unstable, leading to oscillatory behavior. 
A simple numerical trick has been proven very powerful 
in suppressing such oscillations. We introduce a damping 
factor < / < 1 and iterate, starting with = 0: 

p( m ) = V(p(°)+ g ( ro )) (4a) 
q (m+i) = f q (m) /)n p ( m ) (4b) 

We increase / adaptively when oscillations in E(q^) 
are detected and decrease it when convergence is stable. 
All results in Refs. [6,5,9] were obtained using this iter- 
ation procedure, which appears to be both remarkably 
stable and fast. For example, it took just 25 iterations 
to achieve 10 -7 accuracy in the components of q for the 
cluster of 2000 molecules of PTCDA, mentioned above. 

Iterations (4) usually perform quite well. However, 
they are not guaranteed to converge and sometimes fail. 
For example, we found intermittent convergence prob- 
lems for a pentacene cation in a layer on a metal sub- 
strate. Image charges in the metal induce strong charge 
redistribution along the pentacene axis with largest 
polarizability. On the other hand, the surrounding 
molecules in the layer have the opposite effect of drawing 
charge into the inside of the layer, closer to the center of 
the ion. As follows from the self-consistent solution, the 
second effect outperforms the first. This clearly leads to 
difficulties because the initial iterations tend to drive the 
system away from the convergence point, indicating also 
that a perturbative solution may be counterintuitive in 
this case. 

Computational difficulties of this kind are rare. We 
now describe an alternative approach to solving polariza- 
tion equations, which is formally stable and is guaranteed 
to converge, while it may not be as efficient as Eqs. (4). 

B. Subspace methods for polarization problems 

Let us return to the original Eq. (2). We can treat it 
as a general linear system 

Aq = 6, (5) 

with a symmetric matrix A = G + V on the left-hand 
side. Solution of such a linear system is equivalent to 
minimizing the energy functional 

E{q) = \{qAq)-{bq), (6) 



whose minimum exists when A is positive-definite. 

We present a variationally-stable subspace approach to 
solving (5) based on the minimum principle (6). By re- 
scaling the solution we can assume the right-hand side to 
be normalized, (bb) = 1. We shall be searching for the 
solution q in the Krylov subspace K. of the matrix A, gen- 
erated by the vector b: K = span(6, Ab, A 2 b, ...,yl m ~ 1 &). 
At every step m we add a new vector to the subspace 
JC and construct the best vector q which minimizes (6) 
within /C. This guarantees that each step m yields a 
better solution, and thus the procedure converges. 

We use the Hcrmitian Lanczos recursion [10,11] to 
build an orthonormal basis in /C: q^ = b, 

<l {m+1) = ^!+i(M {m) - a m q (m) - P m q {m - iy ) (7) 

where a m and j3 m are chosen at each step m to orthogo- 
nalize g( m+1 ) with respect to both q( m ^ and q^" 1 " 1 ^, en- 
suring that (g("*> \q^) = 5 mn for every m and n. Ex- 
pressing q = ^2 m c m q^ m ^ and solving for the coefficients 
c rn that minimize E(q) we find that they obey a linear 
system A mn c„ = b m , where A mn = (qW AqW) and 
b m = (q^b). 

The algebraic properties of Lanczos recursion (7) yield 
readily all the components of A and b. The matrix 
A is symmetric tridiagonal with the diagonal elements 
A mm — cs-m and sub-diagonal elements A TOim _i = [3 m . 
Due to the ^special choice of q^ = b, all b m = 0, except 
bi = 1 and b<i = ol\. 

The method described above is variationally stable, 
since it is based on a minimum principle. Its drawback, 
however, is that it is generally slower than iterations (4). 
The reason is that Eqs. (4) follow from the dual form of 
Eqs. (3), which is based in turn on a physical insight of 
solving each molecule separately in the external field of 
other molecules. The subspace approach is not based on 
such an insight. In practice, it takes 200 — 300 iterations, 
rather then 25-35, to achieve convergence with this ap- 
proach. Thus we use the subspace approach only when 
iterations fail. 



C. Subspace method for non-positive definite cases 

The subspace method of the above section being ap- 
plied to Eq. (2) deals with the stiffness matrix G. It is 
preferable numerically to work with n, rather than with 
its inverse. The reason is that whenever there is a small 
polarizability, e.g. in the direction normal to the plane of 
a 7r-conjugated molecule, it translates into a small eigen- 
value of II and, correspondingly, into a large eigenvalue 
of G. Numerically, it is much easier to handle near-zero 
than near-infinity. 

In a more drastic situation, G may have negative eigen- 
values, yet the polarization self-consistent equations (3) 
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continue to make sense. For example, we have intro- 
duced atomic polarizabilities a [6] as the difference a— a c 
between the best available molecular a and a c based 
on semiempirical atomic charges. Usually a c underesti- 
mates a by ~ 10%. But the charge-induced polarizability 
occasionally exceeds a in some directions (see e.g. Ta- 
ble I), and the correction a has one or more negative 
eigenvalues. The self-consistent polarization equations 
(3) then give a saddle point of the energy functional, 
which formally has no minimum. 

In order to use the subspace method to deal with II 
rather than G, we rewrite Eqs. (3) in a matrix form in 
such a way that the matrix on the left-hand side is sym- 
metric: 



-V 1 

i n 



9 (0) 





(8) 



This is again a linear system of type (5), but of twice the 
size. The matrix is, however, not positive-definite, and 
so the above method based on the minimization of (6) 
cannot be applied. 

We present a modified subspace approach that is ap- 
plicable to linear systems of type (5) with symmetric ma- 
trix A that is not necessarily positive definite. Instead of 
searching for the minimum of (6) we minimize the resid- 
ual norm. 



Z{q) = \(Aq-b\Aq-b), 



(9) 



which is also equal to half the norm of the gradient, 
dE{q)/dq = Aq — b. Again, expressing the solution as 
1 = Em c ™5'™'' but solving now for c m that minimize 
Z{q) we find another linear system, ^2 n Z mn c n = w m , 
where Z mn = (q^A 2 q^) and w m = (q^Ab).^ 

After some algebra, we find that the matrix Z is sym- 
metric 5-diagonal with the following nonzero elements: 

Z'ram ftm ftm—li ^m,m— 1 — ip^-m "T" O^m— 1 )/^m— 1 ) 

and Z m>m -2 — /3m-iPm-2, while the right-hand side has 
all components zero except w\ = a\ and w 2 = The 
above expressions for the elements of Z mn can be found 
by consequently multiplying Eq. (7) by itself with m re- 
placed with m — 1, m — 2, etc. We notice by inspection 
that Z = A 2 , except for the last row and column. 

Throughout this section we have assumed that all vec- 
tors are projected onto the subspace Q to preserve charge 
conservation relations. We did not discuss precondition- 
ing techniques, such as rescaling the variables, which may 
improve convergence. 



III. POLARIZATION IN PENTACENE CRYSTALS 

A. Gas-phase properties 

We performed calculations of the electronic structure 
of pentacene in gas phase using B3LYP hybrid density 



functional theory with the extensive 6-311++G(d,p) ba- 
sis set available in Gaussian 98 [12]. D2h symmetry was 
assumed and the geometry was separately optimized for 
the neutral molecule, cation, and anion. Unrestricted 
Hartree-Fock (UHF) was used for the ions. Limited spin 
contamination was observed as estimated by the maxi- 
mum deviation of the total spin S 2 = 0.769 before anni- 
hilation. Comparison of the energies of the ions in the 
neutral and ionic optimized geometries yields the relax- 
ation energies A+ = 46 mcV for the cation and A_ = 68 
meV for the anion, under the assumption of equal zero- 
point contributions for the neutral molecule and ions. A + 
compares well with the value of 59 meV derived from ex- 
periment [13]. 

The vertical ionization potential is found to be I = 
6.229 eV, which is significantly below the recently re- 
ported experimental value of 6.589 eV [13]. Other values 
of 6.64 eV [14] and 6.74 eV [15] were reported previously. 
The calculated value of / must be corrected by the differ- 
ence in zero-point energies before the direct comparison 
can be made, but zero-point differences are probably less 
than 360 meV. 

The experimental value for the electronic affinity is 
known with an error bar. Our calculated B3LYP value 
A = 1.475 eV is 105 meV above the recommended av- 
erage 1.37 eV [16] of experimental data. Thus, B3LYP 
values lead to the gas-phase pentacene charge gap I—A = 
4.754 eV (without zero-point correction) vs. 5.22 eV de- 
rived from experiment. 



B. Polarization energies 

Following the procedure in [6] we first calculate 
atom-atom polarizabilities liij for individual pentacene 
molecules using the semiempirical INDO/S Hamiltonian 
[17]. Experimental geometries from the crystal structure 
data [18] were used. Positions of hydrogen atoms, not 
known experimentally, were AMl-optimized using Gaus- 
sian 98 with the heavy atoms fixed. 

Due to the inequivalence in crystal field for the two 
molecules in the unit cell, their geometries are slightly 
different. There are two atom- atom polarizability ten- 
sors Tlij with i, j ranging over 36 atoms of a pentacene, 



IT 



dpi 



d 2 E 



(10) 



E is the INDO/S ground-state energy, fa = <j>{ri) is the 
crystal potential at atom i and pj are Lowdin atomic 
charges. Charge redistribution is then described explic- 
itly as 



j 



(11a) 
(lib) 
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[cf. Eqs (3)], where p\ are gas- phase charges, /J,i are 
induced atomic dipoles and Ft = — V0(rj) are electric 
fields, a — a — a c is a correction [6] that accounts for 
the difference in ab-initio molecular polarizability and the 
"charges-only" part a c associated with : 

ij 

We computed a for the two geometries using B3LYP den- 
sity functional with 6-311++G(d,p) basis set of Gaussian 
98. 

Principal values of a and "charges-only" a c are com- 
pared in Table I. Charge redistribution accounts quite 
well for the in-plane components «ll and auM of the 
neutral molecule. The polarizability qnn, normal to the 
molecular plane, is completely "atomic," as expected, 
and cannot be described by charge redistribution [9,6]. 
Values that are not exactly zero are due to non-planar 
molecules in the crystal. The INDO/S results for ions are 
far too large; we return to this point in the Discussion. 

Table I. Molecular polarizability of pentacene along prin- 
cipal axes (L = long, M — medium, N = normal) 



Molecule cunn (A 3 ) a M M (A 3 ) a LL (A 3 ) 



B3LYP/6-311++G** 



Neutral(l) 


17.77 


37.66 


91.43 


Neutral(2) 


18.02 


38.02 


99.48 


Anion(l) 


19.87 


40.16 


135.47 


Anion(2) 


20.20 


40.81 


131.88 


Cation(l) 


16.31 


35.91 


130.87 


Cation(2) 


16.54 


36.54 


124.14 


a c , Eq. 


(12) 






Neutral(l) 


0.06 


39.26 


87.41 


Neutral(2) 


0.05 


39.61 


95.83 


Anion(l) 


0.06 


35.26 


194.35 


Anion(2) 


0.05 


39.98 


231.24 


Cation(l) 


0.06 


36.03 


368.28 


Cation(2) 


0.05 


40.55 


282.52 



The difference between the B3LYP and a c in Table I is 
distributed over the atoms with the weight proportional 
to the atomic valence charge [6] . Since charge redistribu- 
tion overestimates the in-plane molecular polarizability 
in pentacene, the atomic correction tensor a = a — a c is 
not positive definite. It is still small, however, for neutral 
molecules, and the self-consistent solution is well defined, 
even though it does not correspond strictly to the mini- 
mum of the total energy (cf. discussion in II. C). 



C. Properties of the neutral lattice 

Both pentacene [18] and anthracene [19] illustrate 
herring-bone packing with two molecules per unit cell, 
but their Bravais lattices are different. Pentacene is tri- 
clinic [18]. The two molecules are inequivalent and sub- 
ject to different crystal-field environments. The polariza- 
tion energy of a cation or anion consequently depends on 
which molecular site is charged. Since crystalline electric 
fields are perturbations, we expect comparable polariza- 
tion energies that differ at most by ~ 100 meV. Such dif- 
ferences have already been computed [20] in pentacene 
crystals using the submolecular method [7,21] in which 
the gas-phase molecular polarizability a is placed as a/5 
at the ring centers. 

Solving self-consistent polarization equations for the 
neutral lattice of pentacene, we find that polarization 
contribution to sublimation energy is negligible (1.23 
meV per molecule) , as in anthracene, since the gas-phase 
charges are small due to an approximate electron-hole 
symmetry of the valence shell. We find the dielectric 
tensor of pentacene to be highly anisotropic, with prin- 
cipal values Ki = 5.336, K2 = 3.211, and K3 = 2.413. 
For reference purposes, the directional cosines of the 
principal axes are ni = (—0.296,-0.314,0.902), 112 = 
(-0.021,0.946,0.322), and n 3 = (0.955,-0.076,0.286), 
respectively, in the Cartesian coordinate system with 
the pentacene lattice vectors (in A) a = (7.900,0,0), 
b = (0.444,6.044,0), and c = (-6.153,-2.858,14.502). 
The direction of ni coincides with the direction of long 
axes of pentacene molecules to the accuracy to which the 
latter can be defined (<~ 5°). 

The calculated principal axes and values of k for an- 
thracene crystals agree with experiment [9]. Measure- 
ments of refractive indices or optical dielectric constants 
are challenging in anisotropic molecular crystals [22] . We 
are not aware of such data for pentacene crystals. Di- 
electric data will be necessary for organic devices with 
improved performance. 

D. Charge carriers in the bulk 

We compute the polarization energy of charge carri- 
ers following [6], by placing an ion in the infinite neu- 
tral lattice and considering an imaginary sphere of M 
molecules surrounding it. We allow only the molecules 
whose centers fall into the sphere to relax their charge 
distributions from the self-consistent values for the neu- 
tral lattice and monitor convergence as M is increased. 
We use the same molecular geometry for the ion. The 
small (~ 100 meV) relaxation energies of large aromatic 
molecules make the same geometry both an excellent and 
convenient approximation. Atom-atom polarizabilitics 
Uij for the ion are computed using unrestricted Hartree- 
Fock (UHF) INDO/S. 
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Fig. 2 Convergence of P± with increasing cluster size M. 

Figure 2 shows linear convergence of polarization en- 
ergy with Af" 1 / 3 , which is proportional to the inverse 
radius of the sphere. The polarization energy extrap- 
olates as oc M -1 / 3 , because the missing part at large 
finite M can be thought of as polarization of an infinite 
homogeneous dielectric with a spherical cavity of radius 
R oc M -1 / 3 , -(1 - K~ 1 )e 2 /R [6]. Extrapolation to in- 
finite radius gives the polarization energy in the bulk, 
P = P + + P_ = 2.014 cV and 2.007 eV for the ion oc- 
cupying positions 1 and 2 in the unit cell, respectively. 
The "charges only" values of P = P + + P_ are - 10% 
less, as also found in other organic crystals. Since n e s for 
an anisotropic dielectric medium is known [23] , the slope 
in Fig. 2 at large M is directly related to the n a found 
above, and the slopes agree within a few percent. We 
note that n~g implies identical slopes at large M for the 
anion and cation in either position, as found separately 
in Fig. 2. 

Using B3LYP values for / and A in the gas phase, 
we calculate the transport gap of the pentacene crystal 
E t — I — A — P = 2.740 cV in the limit of zero over- 
lap. This correlates quite well with the reported band 
gap of 2.85 eV, obtained from careful fit of electroab- 
sorption spectra [20]. E t = 2.78 eV was reported in [24] 
also based on electroabsorption spectra of charge-transfer 
(CT) states. The experimental value for the charge gap, 
I — A = 5.22 eV, used in this context results in E t = 3.21 
eV, which is too high. 

E. Ion pairs and CT states 

We next report polarization energies P pa i r of CT states 
with an ion at the origin and the counter ion at a nearby 
site. We place the ion pair within an imaginary sphere 
in a neutral lattice and calculate the effective interaction 
V^ff, which is the difference of the total polarization en- 
ergy and the energy of a pair of well-separated charges, 



V cS = (P pair - P+ - P-). -V eS is the binding energy of 
CT states in the limit of no overlap and is closely related 
to the binding energies of Frenkel excitons, when the final 
state is a molecular excitation. 

Fig. 3 shows convergence with the size of the sphere. 
Since the total ion charge within the sphere is zero, the 
leading term cx M~ 1//3 vanishes. The polarization en- 
ergy of a dipolc in a dielectric medium converges much 
faster as oc M _1 . Extrapolating to infinite M, we ob- 
tain the lowest CT states at -0.719 and -0.679 eV for 
the nearest-neighbors (ccnter-to-center distance 4.799 A) 
and —0.685 and —0.675 cV for next-nearest neighbors 
(5.151 A), as listed in Table II. These results compare 
reasonably with submolecular calculations [20] of —0.777 
and —0.698 eV, respectively, with a/5 at ring centers. 
Submolecular results do not distinguish between anions 
and cations and also depend on the precise partitioning 
of a, which is left open. 
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Fig. 3 Convergence of the effective electron-hole interac- 
tion energy with cluster size M. 

Table II. Effective interaction energies V e fi(r), in eV, 
of cation-anion pairs. The submolecular data is from [20]. 
Center-to-center distances r are also listed for reference. 



Anion Cation 


r, A 


Keff(r) 


T/submol/„\ 
V cS \ r ) 


(0,0,0) (-i,i,0) 
(-i,i,0) (0,0,0) 
(0,0,0) (I, |,0) 
ihho) (0,0,0) 


4.7990 
4.7990 
5.1514 
5.1514 


-0.719 
-0.679 
-0.685 
-0.675 


-0.756 
-0.756 
-0.698 
-0.698 



It is interesting that the mutual orientation of the 
cation and ion together with variations in their electro- 
static response are enough to reverse the order of states, 
such that the pair with larger center-to-center distance 
has greater binding energy (lines 2 and 3 in Table II). 
This is an effect very similar in nature to the one sug- 
gested by Mazur and Petelenz [25]. 



6 



IV. POLARIZATION IN PENTACENE FILMS 

A. Inert metal-organic interfaces 

We now address polarization energy in thin pentacene 
films deposited on metal surface (Fig. 1). We use an ap- 
proach similar to the one sketched recently for thin films 
of PTCDA [5], which we present in greater detail. We 
assume an ideal metal-organic interface without chemi- 
cal interactions. Organic molecules are physisorbed and 
in van der Waals (vdW) contact with the metal. We 
model the metal as an equipotential surface parallel to 
the organic layer. The potential at the metal is assumed 
to be <j) = 0, as the actual value drops out from the 
combined polarization energy for the electron and hole, 
P = P + + P_, even though we perform separate calcu- 
lations for P + and P_. In fact, any additive constant 
in the potential cancels out for any neutral entity with 
an arbitrary but equal numbers of positive and negative 
charges. 

For the organic film we use the bulk crystal geometry, 
assuming a parallel slab of the crystal cut along molecu- 
lar layers and placed parallel to the metal surface. The 
only parameter in the model that allows limited adjust- 
ment is the actual distance of the equipotential plane 
from the slab. We choose it such that the equipotential 
surface is at the vdW contact distance of 2.80 A between 
the closest hydrogen and gold atoms. The molecular ar- 
rangement depicted in Fig. 1 is in fact drawn to scale. 

As in [5] the equipotential metallic surface is treated 
by means of the image charges. For each atomic partial 
charge and induced dipole we assign and place an oppo- 
site charge and dipole in mirror positions when comput- 
ing potentials and potential gradients. The charges and 
dipoles within a molecule do not contribute to the poten- 
tials and potential gradients that act on that molecule. 
However, the images of these charges and dipoles do con- 
tribute. For example, a single ion exhibits no polarization 
energy in gas phase, while it acquires polarization energy 
when placed near the metal surface. 

As in bulk calculations, we solve the polarization prob- 
lem for slabs in two steps. We first consider a neutral 
film of N molecular layers with no ions and with trans- 
lational symmetry in two dimensions. We solve the self- 
consistent equations and find ground-state charges 
and dipoles Jt^ , which are different for each molecular 
layer k — 1, N and for each type j = 1, 2 of molecule 
in the unit cell. While the slab is infinite in two dimen- 
sions, the number of variables is finite due to translational 
symmetry. 

We then replace one molecule in the surface molecular 
layer by an ion and solve for the difference Spf = pi —~pj 
and 6/j,f = fif — Ji i 3 for every molecule a in the slab. 
Since the slab is infinite, we define clusters whose shapes 
are pillboxes with variable radius R parallel to the metal 



and fixed thickness to include N molecular layers and 
their images. We then relax the charge distribution of 
all M molecules whose centers are in the pillbox. We 
assume 6p1 = <5/x" = for the molecules a outside the 
pillbox, which means that their charge distributions are 
not relaxed. They contain, however, partial charges and 
dipoles 'p i ■* and Ji i ■* , which contribute to the total polar- 
ization energy of the film. 

B. Surface and interface polarization 

The polarization energy for the ion appears, as in the 
bulk, as a finite difference of two infinite energies: for the 
infinite slab with and without the ion. The expression for 
the energy that takes care of this cancellation is the same 
as Eq. (27) in [6], but with the potentials and potential 
gradients containing contributions from image charges 
and dipoles. We perform calculations for finite pillbox 
clusters of increasing M, and extrapolate to M — ► oo. 

Fig. 4 shows the convergence of polarization energies 
for cations (upper panel) and anions (lower panel), as 
well as the extrapolated values for pentacene films one 
to ten molecular layers thick. Two lines for each film 
thickness correspond to two types of molecules in the 
unit cell. Convergence with increasing pillbox radius R 
is shown as I /Mi, where M; = M/N is the average num- 
ber of pentacenes per layer. The computational effort 
increases with N. The largest N = 10 clusters in Fig. 3 
contain M — 2806 molecules, each containing 36 atoms, 
and their images. 

The combined extrapolated values P s = Pf + P^ are 
also listed in Table III and plotted vs. 1/N in Fig. 5. 
The monolayer (N — 1) and bulk values are almost the 
same. Image charges and dipoles at vdW separation are 
sufficiently more polarizable to compensate for a vacuum 
on the other side. P S (N) decreases with increasing N 
as less polarizable organic layers are introduced between 
the surface and the metal. In the limit of a thick film, or 
a free pentacene ab surface, we estimate P s — P ~ —0.23 
eV. 

Table III. Variation in polarization energy at the surface 
(P s ) and at the metal/organic interface (P M ) in an N- layer 
thick pentacene film on a metallic substrate. The values are 
reported for anion at position (1) and cation at position (2) 
in the unit cell, which correspond to the lowest energy state. 



Layers, N 


(P s -P), mcV 


(P M - P), meV 


1 


6 


6 


2 


-40 


113 


3 


-71 


125 


5 


-108 


129 


10 


-166 


130 


thick film 


-227 


130 
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In Table III the interface values P M = Pf + P M refer 
to the polarization energy of an ion next to the metal in 
a film of N molecular layers. The biggest increase of 100 
meV at iV = 2 comes from the first pentacene overlayer. 
Additional layers produce small changes and the inter- 
face value in a thick films is P M — P = 0.13 eV. The hole 
and electron components of P M are relevant for matching 
transport levels to Fermi energies for facile injection of 
carriers. Such optimization for organic devices is largely 
empirical at present, often without any consideration of 
polarization. The variation of P s (N) or P M (N) with film 
thickness is less in pentacene than in PTCDA [5]. Pen- 
tacene layers are several times thicker, with the long axis 
normal to the metal, and a single layer already provides 
an effective polarization environment for the charge. 
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Fig. 4 Convergence of polarization energies of a cation 
(upper panel) and anion (lower panel) in the outermost sur- 
face layer of pentacene film on gold vs. the average number 
Mi = M/N of molecules in each molecular layer (M molecules 
in N layers in the pillbox). Open and solid symbols corre- 
spond to the molecules of type (1) and (2) respectively. 
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Fig. 5 Variation in polarization energy at the surface 
(P S (N) -P) and at the metal/organic interface (P M (N) - P) 
for an iV-layer thick pentacene film on a metal substrate. 



C. Molecular potential 

We have already noted the competition between im- 
age charges and induced dipoles in films of molecules 
with large polarizability normal to the surface. The self- 
consistent potential </> a (r) at molecule a is due to the 
polarized densities p b (r) of all other molecules and their 
images. A single pentacene cation in the surface layer in- 
duces an image anion in the metal. Since the interaction 
is attractive, the cation's <f> cat (r) increase monotonically 
with the distance from the surface. Similarly, the <fi an (r) 
of a single anion on the surface decrease monotonically 
with distance from the surface. 

The situation in films is quite different due to the sur- 
rounding neutral molecules. Since the large pentacene 
polarizability is normal to the surface, the induced mo- 
ments due to the image anion are parallel. The repulsion 
of parallel dipoles and image dipoles is relieved by re- 
distributing charge towards the middle of the molecular 
layer. Now the (j) cat (r) has a minimum at atoms near the 
center, while the t/>°"(r) has a maximum near the cen- 
ter. The actual <f) zon (r) in the limit of zero overlap is the 
solutions of the self-consistent equations. 

Figure 6 shows the self-consistent potential </> ca '(r) for 
r along the axis of a cation in a pentacene monolayer 
on metal. The gradients yield electric fields ±10 7 V/cm 
that reverse over a single molecule, directly confirming 
that fields are strong and non-uniform. The variation 
of <fi cat (r) with r shows the importance of repulsion be- 
tween induced dipoles. Image charges account for more 
negative 4> cat (r) close to metal. Similar variations of self- 
consistent potentials are found at ions (1) and at neigh- 
boring neutral molecules. The strong shielding of a single 
layer of pentacene follows from the weak dependence on 
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N. Similar curves (not shown) for N = 2, 3, and 5 are 
essentially parallel to the N = 1 curve in Fig. 6, but 
displaced to less negative (positive) potential for cations 
(anions). 



-0.8 




-2.4 1 ' ' ' 1 

5 10 15 20 

Distance from metal, A 

Fig. 6 Self-consistent potential for a cation created by the 
polarization field of neutral molecules in the lattice and the 
image charge distribution in the metal. Data shown for a 
cation in a monolayer on metal substrate, placed in the po- 
sition (2) in the unit cell. The sketch of the molecule shows 
actual position with respect to the metal. 

V. DISCUSSION 

The bulk polarization energy can be estimated simply 
in terms of molecular size, P± = (e 2 /2R)(l — 1/k), for a 
sphere of radius R in an isotropic dielectric medium with 
k ~ 3 for organic molecular crystals. Such rough values 
rationalize much data but preclude accurate positioning 
of transport or CT states. The effective pentacene ra- 
dius is 14% greater than that of anthracene, while the 
polarization energy decreases by only 9%. On the other 
hand, the radii of pentacene and PTCDA molecules are 
the same within 3%, but the polarization energy differs 
by 10%. The greatest contrast to PTCDA films [5], is 
that the polarization energy P s is within 10% of the 
bulk value in pentacene, compared to P s / P — 0.77 in 
PTCDA. 

Direct relations between gas-phase and crystalline en- 
ergies have long been sought using molecular exciton the- 
ory [7,26]. Accurate evaluation of electronic polarization 
permits sharper comparisons. There is considerable scat- 
ter in reported values for both the gas and solid, as noted 
in Sections III(A) and HI(D) for the ionization potential 
and electron affinity of pentacene. The transport gap is 
I — A — P in the limit of zero intermolccular overlap. 
The experimental values for the pentacene transport gap 
range from 2.85 eV [20] to 2.4 eV suggested in [27], to 
2.2 eV used in [28]. This scatter was not "important" 



when polarization energies had to be estimated in the 
bulk and were not available at all at surfaces. With the 
ability to calculate polarization energies, it becomes pos- 
sible to analyze and connect various sources of experi- 
mental data on such quantities as optical and transport 
gaps, CT binding energies, and work functions. 

All molecular calculations, except for the gas-phase 
properties, were performed using experimental geome- 
tries from the crystal structure data [18]. Slight varia- 
tions in the geometry of the two molecules in the unit 
cell were consequently taken into account. While the 
variations in electronic energies are small (tens of meV), 
which is consistent with the adopted approximation, we 
believe that crystalline geometries are preferable to gas- 
phase geometries when vertical electronic transitions are 
considered in the condensed phase. 

We note that INDO/S greatly overestimates the in- 
plane polarizabilities of ions in Table I. The huge differ- 
ence of a between cations 1 and 2 is clearly unphysical 
in view of small distortions in the crystal, and it rational- 
izes why convergence was particularly difficult for cation 
1. The UHF INDO/S approach fails for ions. Radical 
ions of large conjugated molecules are expected to be 
more polarizable than the neutral molecules, in line with 
the B3LYP entries in Table I. While the idea is that a 
constitutes a small correction, and it is not so for ions, 
we expect this to have a negligible impact on the results. 
The reason is that the largest part of the polarization en- 
ergy comes not from the ion itself but from the polariza- 
tion of many neutral molecules surrounding it. Moreover, 
the molecular a is kept no matter how it is partitioned 
between charge redistribution (a c ) and induced dipoles. 
In the future we may choose to use Hy of the neutral 
molecule for the ions, and to absorb the polarizability 
difference in a. This requires separate B3LYP calcula- 
tion of a for the ions, as done here, but not in [6] . 

Small atomic charges in acenes are related to approx- 
imate electron-hole symmetry. The crystal potential 
4>i = 4>( r i) at atom i due to INDO/S charges is con- 
sequently small, as shown by ~ 1 meV per molecule con- 
tributions to the sublimation energy. Since the potential 
due to the B3LYP or any other gas-phase charge distri- 
bution can readily be computed at the positions of other 
atoms in the crystal, we can obtain first-order correc- 
tions to polarization energies in general [29] by combining 
the best available gas-phase potential with self-consistent 
charges and dipoles based in INDO/S calculations. The 
quadrupole moments of acenes lead to opposite shifts 
of P + and P_, as discussed previously for submolccules 
[7,21]. Since such corrections cancel in P = P + + P_ 
in anthracene [6], we expect them to be small in pen- 
tacene. They are not negligible in PTCDA or conju- 
gated molecules that contain heteroatoms and hence par- 
tial charges. 

Our values for the the transport gap are likely an over- 
estimate, due to the bandwidth effects for electrons and 
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holes, and the correction can be introduced as E t — 
I — A — P — W . The bandwidth W in molecular crystals 
at room temperature is thought to be of the order of few 
tens of meV or less, which justifies our zero-overlap ap- 
proximation. At lower temperatures, however, the band- 
widths may be large in certain directions. It has been 
suggested that the bandwidth in pcntacene (at zero tem- 
perature) may be anomalously large and reach 500 meV 
[30]. If we estimate W ~ 0.5 eV, we get E t = 2.71 eV 
with the experimental charge gap, while the calculated 
charge gap leads to the value 2.35 eV, which is too small. 

When bandwidths arc comparable to polarization en- 
ergies, the zero-overlap approximation cannot be justi- 
fied, and we face a challenging problem of combining the 
"localized" treatment of polarization effects with the "de- 
localized" language of band theory. The self-consistent 
polarization approach can provide a proper "zero-order" 
approximation for the attempts to treat bandwidth ef- 
fects as perturbation. 

We note that our self-consistent approach has the abil- 
ity to distinguish between the electronic polarization of 
cations and anions, since the quantum-mechanical struc- 
ture of molecular response is used to obtain the redis- 
tribution of charge. We can also estimate variations of 
polarization caused by slight changes in molecular geom- 
etry between crystal and gas-phase, and especially be- 
tween inequivalent molecules in the unit cell. These vari- 
ations appear to be of the order of 20 — 40 meV, which is 
on the threshold of the accuracy, but such variations are 
probably correct on the order of magnitude, and may be 
important in the analysis of experimental data. 

In summary, we have applied the recently developed 
self-consistent approach based on evaluating charge re- 
distribution in organic molecules to the problem of elec- 
tronic polarization in pentacene molecular crystals. The 
power of the approach comes from the combination of 
semiempirical treatment of charge distribution by means 
of atom-atom polarizability tensors with accurate ab ini- 
tio gas- phase calculations of molecular polarizabilitics. 
As a result, we are now able to calculate polarization 
energies with accuracy on the order of 0.1 eV or better. 
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